Energy cannot be created or destroyed, therefore it must be
conserved. To see this, consider the following system defined as a
volume $V$ in space enclosed by the surface $S$ as shown in figure
\ref{fig:ConsEnergy}. The system, is filled with a material
characterized by the constitutive parameters $\varepsilon$, $\mu$,
$\sigma_e$ and $\sigma_m$. Inside the region are the current
sources $\vec{\mathcal{J}_i}$ and $\vec{\mathcal{M}_i}$ that
produce the electric field $\vec{\mathcal{E}}$ and the magnetic
field $\vec{\mathcal{H}}$. These fields must obey Maxwell's
equations.

%\includegraphics[angle=0, width=1.0\textwidth]{./Figures/cavity_junction.eps}
\begin{figure}[!hbp]
\begin{center}
\includegraphics[width=0.5\textwidth,height=0.5\textwidth]{./Figures/ConsEnergy.eps}
\end{center}
\caption{The system, defined as a volume of space bounded by the
surface S is filled with a material characterized by the
constituent parameters
($\sigma_m$,$\sigma_e$,$\mu$,$\varepsilon)$. Within the system are
current sources [$\vec{\mathcal{J}_i}$,$\vec{\mathcal{M}_i}$] that
produce the fields
[$\vec{\mathcal{E}}$,$\vec{\mathcal{H}}$,$\vec{\mathcal{D}}$,$\vec{\mathcal{B}}$].}
\end{figure}\label{fig:ConsEnergy}

In order to derive the \emph{conservation of energy} for
electromagnetics, equation (\ref{eqn:fl}) is scalar multiplied by
$\vec{\mathcal{H}}$ and (\ref{eqn:al}) is scalar multiplied by
$\vec{\mathcal{E}}$,
\begin{align}
    \vec{\mathcal{H}}\cdot\nabla\times\;\vec{\mathcal{E}}&=-\vec{\mathcal{H}}\cdot\left(\frac{\partial\vec{\mathcal{B}}}{\partial{t}}+\vec{\mathcal{M}}\right)\label{eqn:fldot}\\
    \vec{\mathcal{E}}\cdot\nabla\times\,\vec{\mathcal{H}}&=\vec{\mathcal{E}}\cdot\left(\frac{\partial\vec{\mathcal{D}}}{\partial{t}}+\vec{\mathcal{J}}\right)\label{eqn:aldot}
\end{align}
Subtracting (\ref{eqn:fldot}) from (\ref{eqn:aldot}) yields,
\begin{equation}\label{eqn:consofenergy1}
    \vec{\mathcal{H}}\cdot\nabla\times\;\vec{\mathcal{E}}-\vec{\mathcal{E}}\cdot\nabla\times\,\vec{\mathcal{H}}=-\vec{\mathcal{H}}\cdot\left(\frac{\partial\vec{\mathcal{B}}}{\partial{t}}+\vec{\mathcal{M}}\right)-\vec{\mathcal{E}}\cdot\left(\frac{\partial\vec{\mathcal{D}}}{\partial{t}}+\vec{\mathcal{J}}\right)
\end{equation}
Using the vector identity,
$\nabla\cdot(\vec{A}\times\vec{B})=\vec{B}\cdot\nabla\times\vec{A}-\vec{A}\cdot\nabla\times\vec{B}$
reduces (\ref{eqn:consofenergy1}) too,
\begin{equation}\label{eqn:consofenergydiff2}
    \nabla\cdot(\vec{\mathcal{E}}\times\vec{\mathcal{H}})=-\vec{\mathcal{H}}\cdot\left(\frac{\partial\vec{\mathcal{B}}}{\partial{t}}+\vec{\mathcal{M}}\right)-\vec{\mathcal{E}}\cdot\left(\frac{\partial\vec{\mathcal{D}}}{\partial{t}}+\vec{\mathcal{J}}\right)
\end{equation}
Substituting in (\ref{eqn:Js}) and (\ref{eqn:Ms}) and rearranging
yields,
\begin{equation}\label{eqn:consofenergydiff}
    \nabla\cdot(\vec{\mathcal{E}}\times\vec{\mathcal{H}}) + \left(\vec{\mathcal{H}}\cdot\frac{\partial\vec{\mathcal{B}}}{\partial{t}}+\vec{\mathcal{E}}\cdot\frac{\partial\vec{\mathcal{D}}}{\partial{t}}\right) + \left(\vec{\mathcal{H}}\cdot\vec{\mathcal{M}_i}+\vec{\mathcal{E}}\cdot\vec{\mathcal{J}_i}\right) + \left(\vec{\mathcal{H}}\cdot\vec{\mathcal{M}_\ell}+\vec{\mathcal{E}}\cdot\vec{\mathcal{J}_\ell}\right)=0
\end{equation}
Using the \emph{divergence theorem} of
(\ref{eqn:divergencetheorem}), equation
(\ref{eqn:consofenergydiff}) can be changed to integral form,
\begin{multline}\label{eqn:consofenergyint}
    \oint_S\vec{\mathcal{E}}\times\vec{\mathcal{H}}\cdot{d\vec{s}} +
    \left(\int_V\vec{\mathcal{H}}\cdot\frac{\partial\vec{\mathcal{B}}}{\partial{t}}\;dv+\int_V\vec{\mathcal{E}}\cdot\frac{\partial\vec{\mathcal{D}}}{\partial{t}}\;dv\right)\\ +
    \left(\int_V\vec{\mathcal{H}}\cdot\vec{\mathcal{M}_i}\;dv+\int_V\vec{\mathcal{E}}\cdot\vec{\mathcal{J}_i}\;dv\right) +
    \left(\int_V\vec{\mathcal{H}}\cdot\vec{\mathcal{M}_\ell}\;dv+\int_V\vec{\mathcal{E}}\cdot\vec{\mathcal{J}_\ell}\;dv\right)
    =0
\end{multline}
Equations (\ref{eqn:consofenergydiff}) and
(\ref{eqn:consofenergyint}) are interpreted respectively as the
\emph{conservation of energy} equations in differential and
integral form. The quantity
$\vec{\mathcal{E}}\times\vec{\mathcal{H}}$, known as the
\emph{Poynting vector}, is recognized as a power density and has
units of $\frac{V}{m}\frac{A}{m}=\frac{W}{m^2}$. From
(\ref{eqn:consofenergyint}), the total power exiting the volume
$V$ through the surface $S$ is,
\begin{equation}\label{eqn:powerexiting}
    \mathcal{P}_e=\oint_S\vec{\mathcal{E}}\times\vec{\mathcal{H}}\cdot{d\vec{s}}
\end{equation}
where $d\vec{s}=\hat{n}\;ds$ and $\hat{n}$ is a unit normal vector
to the surface $S$ directed out from the closed surface. The other
terms in (\ref{eqn:consofenergyint}) have units of power and
energy. The following term from (\ref{eqn:consofenergyint}),
\begin{align}
    \mathcal{P}_s&=-\int_V\vec{\mathcal{H}}\cdot\vec{\mathcal{M}_i}\;dv-\int_V\vec{\mathcal{E}}\cdot\vec{\mathcal{J}_i}\;dv\notag\\
    &=-\int_V\vec{\mathcal{H}}\cdot\vec{\mathcal{M}_i}+\vec{\mathcal{E}}\cdot\vec{\mathcal{J}_i}\;dv\label{eqn:powersupplied}
\end{align}
represents the total power supplied to the system by the impressed
current sources $\vec{\mathcal{M}_i}$ and $\vec{\mathcal{J}_i}$.
The term,
\begin{align}
    \mathcal{P}_d=\int_V\vec{\mathcal{H}}\cdot\vec{\mathcal{M}_\ell}\;dv+\int_V\vec{\mathcal{E}}\cdot\vec{\mathcal{J}_\ell}\;dv\label{eqn:powerdissipated0}
\end{align}
is the total power dissipated from the system by the lossy current
densities $\vec{\mathcal{M}_\ell}$ and $\vec{\mathcal{J}_\ell}$.
These are sometimes referred to as the system sinks (opposite of
sources). Substituting (\ref{eqn:MH}) and (\ref{eqn:DE}) and using
the convolution property $a(f\ast g)=(af)\ast g=f\ast(ag)$
%\cite{Bracewell:1999},
\begin{align}
    \mathcal{P}_d&=\int_V\vec{\mathcal{H}}\cdot(\sigma_m\ast\vec{\mathcal{H}})\;dv+\int_V\vec{\mathcal{E}}\cdot(\sigma_e\ast\vec{\mathcal{E}})\;dv\notag\\
    &=\int_V\sigma_m\ast(\vec{\mathcal{H}}\cdot\vec{\mathcal{H}})\;dv+\int_V\sigma_e\ast(\vec{\mathcal{E}}\cdot\vec{\mathcal{E}})\;dv\notag\\
    &=\int_V\sigma_m\ast\mathcal{H}^2+\sigma_e\ast\mathcal{E}^2\;dv\label{eqn:powerdissipated}
\end{align}
where $\mathcal{E}$ and $\mathcal{H}$ are the magnitudes of the
vectors $\vec{\mathcal{E}}$ and $\vec{\mathcal{H}}$, respectively.
Equation (\ref{eqn:powerdissipated}) represents the total power
dissipated in the system. Energy being stored in the system can be
derived from the following term in (\ref{eqn:consofenergyint}):

\begin{equation}\label{eqn:Pw1} \mathcal{P}_w=\int_V\vec{\mathcal{H}}\cdot\frac{\partial\vec{\mathcal{B}}}{\partial{t}}\;dv+\int_V\vec{\mathcal{E}}\cdot\frac{\partial\vec{\mathcal{D}}}{\partial{t}}\;dv
\end{equation}

Substituting in equations (\ref{eqn:BH}) and (\ref{eqn:DE})
yields,

\begin{equation}\label{eqn:Pw2}    \mathcal{P}_w=\int_V\vec{\mathcal{H}}\cdot\frac{\partial}{\partial{t}}(\mu\ast\vec{\mathcal{H}})\;dv+\int_V\vec{\mathcal{E}}\cdot\frac{\partial}{\partial{t}}(\varepsilon\ast\vec{\mathcal{E}})\;dv
\end{equation}

Using the convolution property of
$\frac{\partial}{\partial{t}}({f}\ast{g})=\frac{\partial{f}}{\partial{t}}\ast{g}=f\ast\frac{\partial{g}}{\partial{t}}$
on (\ref{eqn:Pw2}) yields,

\begin{equation}\label{eqn:Pw3}   \mathcal{P}_w=\int_V\vec{\mathcal{H}}\cdot\left(\mu\ast\frac{\partial\vec{\mathcal{H}}}{\partial{t}}\right)\;dv+\int_V\vec{\mathcal{E}}\cdot\left(\varepsilon\ast\frac{\partial\vec{\mathcal{E}}}{\partial{t}}\right)\;dv
\end{equation}

Applying the convolution property of $a(f\ast g)=(af)\ast
g=f\ast(ag)$ and rearranging yields,

\begin{align}    \mathcal{P}_w&=\int_V\mu\ast\left(\vec{\mathcal{H}}\cdot\frac{\partial\vec{\mathcal{H}}}{\partial{t}}\right)\;dv+\int_V\varepsilon\ast\left(\vec{\mathcal{E}}\cdot\frac{\partial\vec{\mathcal{E}}}{\partial{t}}\right)\;dv\notag\\   &=\int_V\mu\ast\frac{1}{2}\frac{\partial\mathcal{H}^2}{\partial{t}}\;dv+\int_V\varepsilon\ast\frac{1}{2}\frac{\partial\mathcal{E}^2}{\partial{t}}\;dv\notag\\
    &=\frac{\partial}{\partial{t}}\int_V\frac{1}{2}\mu\ast{H}^2\;dv+\frac{\partial}{\partial{t}}\int_V\frac{1}{2}\varepsilon\ast{E}^2\;dv\notag\\
    &=\frac{\partial}{\partial{t}}(\mathcal{W}_m+\mathcal{W}_e)\label{eqn:energystored}
\end{align}

where $\mathcal{W}_m$ and $\mathcal{W}_e$ are respectively the
magnetic and electric energies stored in the system.

The integrand in (\ref{eqn:powerexiting}) represents power density
with units of $(\frac{W}{m^2})$, which can be written as,
\begin{align}
{p}_e&=\vec{\mathcal{E}}\times\vec{\mathcal{H}}\cdot\hat{n}\\
\end{align}
The integrands in (\ref{eqn:powersupplied}) and
(\ref{eqn:powerdissipated}) represents power density with units of
$(\frac{W}{m^3})$, which can be written as,
\begin{align}
p_s&=-(\vec{\mathcal{H}}\cdot\vec{\mathcal{M}_i}+\vec{\mathcal{E}}\cdot\vec{\mathcal{J}_i})\\
p_d&=\sigma_m\ast\mathcal{H}^2+\sigma_e\ast\mathcal{E}^2\\
\end{align}
The electric and magnetic stored energy densities of
(\ref{eqn:energystored}) with units of $(\frac{J}{m^3})$ are
defined as,
\begin{align}
w_m&=\frac{1}{2}\mu\ast{H}^2\\
w_m&=\frac{1}{2}\varepsilon\ast{E}^2
\end{align}
Using the above definitions, the \emph{conservation of energy}
equation (\ref{eqn:consofenergyint}) can be written as,
\begin{equation}\label{consofenergyint2}
    \mathcal{P}_s=\mathcal{P}_e+\mathcal{P}_d+\frac{\partial}{\partial{t}}(\mathcal{W}_m+\mathcal{W}_e)
\end{equation}
which says that the total power supplied $(\mathcal{P}_s)$ is
equal to the summation of the total power exiting
$(\mathcal{P}_e)$ the system plus the total power dissipated
$(\mathcal{P}_d)$ from the system plus the time rate of change
(positive if increasing) of the stored energies
$(\mathcal{W}_e,\mathcal{W}_m)$ of the system.
